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Abstract: 

We study the scaling behavior of the 4D SU(3) lattice gauge theory in the presence 
of a 9 term, by Monte Carlo simulations computing the topological properties at 
imaginary 9. The numerical results provide a good evidence of scaling in the contin- 
uum limit. The imaginary 9 dependence of the ground-state energy turns out to be 
well described by the first few terms of related expansions around 9 = 0, providing 
accurate estimates of the first few coefficients, up to 0(9 6 ). 
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Four-dimensional SU(N) gauge theories have a nontrivial dependence on the 
parameter 9 which appears in the Euclidean Lagrangian as 

to = \f; w (x)f; w (x) - ie^ eflupa F^(x)F^(x), (1) 

where 

?(*) = g^c w i^(x)i^(x) (2) 

is the topological charge density. The 9 dependence of the ground-state energy- 
density F(9), defined as 



Z = J [dA] exp J d A xC^j = exp[-VF(9)), 



(3) 



where V is the space-time volume, can be written as 



m^m-m = \x0 2 s{9) ) (4) 

where x is the topological susceptibility at 9 = 0, 

x = J d A x(q(x)q(0)) e=0 = ^ y~° , Q = J d 4 xq(x), (5) 

and s(9) is a dimensionless even function of 9 such that s(0) = 1. Assuming analyt- 
icity at 9 = 0, s(9) can be expanded as 

s(9) = 1 + b 2 9 2 + 6 4 ^ 4 + • ■ ■ , (6) 

where only even powers of 9 appear. Large- N scaling arguments |1|, 0, |3|] applied to 
the Lagrangian (|]) of a 4D SU(iV) gauge theory indicate that the relevant scaling 
variable in the large- N limit is 9 = 9/N. This implies that in the large- iV limit 
X = 0(1), while the coefficients 6 2i are suppressed by powers of N, i.e. [Q 

b 2i = 0(N~ 2i ). (7) 

The presence of the 9 term has important phenomenological consequences, since 
it violates parity and time reversal symmetry. Experimental bounds on the 9 parame- 
ter in QCD are best obtained from the electric dipole moment of the neutron |5|, [|] , 
which leads to an unnaturally small value for 9, \9\ < 1CT 10 . This suggests the idea 
that there must be a mechanism responsible for suppressing the value of 9 in the 
context of QCD. However, the issue of the 9 dependence, even within pure gauge 
theory, has phenomenological relevance because it provides an explanation to the 
so-called U(1)a problem |7], §|, 0, i-e. explaining the large mass of the r/' meson with 
respect to the 7r meson. The 9 dependence of 4D SU(iV) gauge theories is particu- 
larly interesting in the large- iV limit where the issue may also be addressed by other 
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approaches, such as AdS/CFT correspondence applied to nonsupersymmetric and 
non conformal theories, see e.g. Refs. |TI| [TTj, [12], [13], |I4| . 

Due to the nonperturbative nature of the 9 dependence, quantitative assessments 
of this issue have largely focused on the lattice formulation of the theory, using Monte 
Carlo (MC) simulations. However, the complex nature of the 9 term in the Euclidean 
QCD Lagrangian prohibits a direct MC simulation at 9 ^ 0. Information on the 9 
dependence of physically relevant quantities, such as the ground state energy and 
the spectrum, has been obtained by computing the coefficients of the corresponding 
expansions around 8 = 0. The coefficients of expansion of s(6), cf. Eq. (f|), can 
be determined from appropriate zero-momentum correlation functions @, |15| of the 
topological charge density at 9 = 0, which are related to the moments of the 8 = 
probability distribution P{Q) of the topological charge Q. Indeed 

b * = -^ X4 = ^[(Q 4 >-3(Q 2 > 2 k , (8) 
h = -3^, Xe = - [<Q 6 > - 15(Q 2 )(Q 4 ) + 30(Q 2 ) 3 ] e=0 , (9) 

etc. They parameterize the deviations of P{Q) from a simple Gaussian behavior. 
As shown in Ref. [|16[], the correlation functions involving multiple zero- momentum 
insertions of the topological charge density can be defined in a nonambiguous, reg- 
ularization independent way, and therefore the expansion coefficients &2i are well 
defined renormalization group invariant quantities. The numerical evidence for a 
nontrivial ^-dependence, obtained through MC simulations of the lattice formula- 
tion, appears quite robust. We refer the reader to Ref. [|J for a recent review. On 
the other hand, MC simulations at 9 = have only made it possible to estimate 
the ground-state energy up to the 0{9 A ) term, because the statistical errors rapidly 
increase with increasing order of the expansion, essentially for importance sampling 
problems. The large- N prediction 62 = 0(N~ 2 ) has been already supported by nu- 



merical results [[15|, [IT], [18]. The calculation of the higher-order terms would provide 
a further check of the large- A arguments, which predict them to be suppressed by 
higher powers of N, as in Eq. (^). 

In this paper we consider imaginary values of 9, which make the Euclidean 
Lagrangian ([!]) real, thus allowing us to investigate the imaginary 9 dependence by 
MC simulations. Assuming analyticity at 9 = 0, the results may provide quantitative 
information on the expansion around 9 = 0, cf. Eqs. (|J) and @. Indeed, fits of 
the data to polynomials of imaginary 9 may provide more accurate estimates of 
the coefficients, overcoming the rapid increasing of the statistical errors observed 
at 9 = 0. Another interesting issue concerns the continuum limit of the theory for 
generic, and in particular imaginary, 9 values. Perturbative renormalization-group 
(RG) arguments [pij indicate that 9 is a RG invariant parameter of the theory, thus 
the continuum limit should be approached while keeping 9 fixed to any complex 
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value. Although reflection positivity [20|, and therefore a well defined relativistic 
field theory, requires real 9 values, imaginary 9 may still give rise to a well defined 
continuum limit in the sense of a statistical field theory, with the same asymptotic 
behavior controlled by the perturbative function which does not depend on 9. As 
we shall see, this is indeed supported by the numerical data for the 4D SU(3) lattice 
gauge theory, at least for \9\ < tt, which are well described by the first few nontrivial 
terms of the expansion around 9 = 0. 

Introducing the real parameter 9{, defined by 

9 = -i9i, (10) 

Eq. (||) becomes 
We thus obtain 

X 9 l (l-2b 2 9* + 3b 4 9t + ---), (12) 



(Qhi _ 




V 


d9i 


(Q%. _ 


(Q 2 )e> ~ 


V 


V 



d9' t 
etc. 



X (1 - 6M- + 15Mf + •-•)> (13) 



The 4D SU(iV) gauge theory with an imaginary 9 term can be nonpertubatively 
formulated on the lattice by 

Z L = J [dA] exp (S L + 9 L Q L ) , (14) 
^ = -Tt E ReTrn^(x), /3 = ^, (15) 

where g is the bare coupling, Yl^ u is the standard plaquette operator, given by the 
product of link variables along a 1 x 1 plaquette of the lattice 

IV(x) = U^{x)U v {x + fi + u)Ul(x + u)Ul(x), (16) 

and the topological term 

Q L = J2q L (x) (17) 



is constructed using the lattice operator [21 



±4 

1 

qL{x 



1 ±4 

-55-2 X) w- Tr [ n ^ n P-]- ( 18 ) 



2 4 x 32vr 2 

Notice that this is not the only possible choice for ; the only requirement is that 
it must have the correct continuum limit (0) when a — > (a is the lattice spacing). 



- 3 - 



Our choice of the operator qL,{x) is also motivated by its relative simplicity in view 
of the MC simulations. In the continuum limit |2lJ being a local operator, 

behaves as 

q L (x) — > a% ? (x) + 0(a 6 ), (19) 

where Z g is a finite function of the bare coupling, going to one in the limit ft — > oo. 
Thus, we have the correspondence 

^ = % (20) 

apart from 0(a 2 ) corrections. The renormalization Z q may be evaluated by MC 
simulation at 9 = 0, computing 

7 _ (QQl}9=0 / on 

^"ToW (21) 



where Q is a topological estimator such as those obtained by the overlap method [23 
21 1 or the cooling method [|25|| , which are not affected by lattice renormalizations, 
and, more importantly, by nonphysical background contributions arising when the 
correlations of topological charge densities are measured at coincident points. Note 



that Eq. ( 21 ) assumes that such contributions are also absent when we consider 
expectation values of mixed products of Q and Ql', this hypothesis will be checked 
by the scaling consistency of the results. In particular, we expect that the ratios 

;<7M - 6i - 2b 2 9f + 36 4 0? + (22) 
1-6M- +15&40* + -, (23) 



(Q 2 )e=o 



have a well defined continuum limit as functions of Q{ . 

We mention that similar approaches based on MC simulations with imaginary 9 
terms have been also pursued to study the 6 dependence of 2D CP^ -1 models, see 



Refs. H, |3,|2J, [2J. 

In the following we present results of MC simulations of the 4D SU(3) lattice 
gauge theory, at (3 = 5.9, 6, 6.2, for lattice sizes L = 16, 16, 20, respectively; the 
simulations are carried out both at 6l = and 9l ^ 0, within the region \8i\ < ir. 
The updating of the link variables is performed using the over relaxation algorithm 



described in Ref . |3(J • Since our numerical study requires high-statistics MC simula- 
tions, we choose the cooling method as estimator of the topological charge Q, and in 
particular the implementation outlined in Ref. fT5 |. The topological charge has been 



measured on cooled configurations (by locally minimizing the lattice action), using 
the twisted double plaquette operator (|T8|). As is well known, the sum over the whole 
lattice of the twisted double plaquette, measured on cooled configurations, takes on 
values Q t — ka, where k is an integer and a < 1. We determine the typical value of 
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L 


x = (Q 2 ) e =o/V 


b 2 


6 4 


z q 


5.9 


16 


0.0001532(4) 


-0.026(6) 


-0.016(13) 


0.1122(6) 


6.0 


16 


0.0000743(2) 


-0.027(4) 


0.004(3) 


0.1353(9) 


6.2 


20 


0.0000201(1) 


-0.028(4) 


0.000(2) 


0.174(2) 



Table 1: 9 = MC results for the 4D SU(3) lattice gauge theory. 






Link to update 
Force 



Figure 1: The force corresponding to a given link. All but the first contribution arise 
from the 6 term. 

a by minimizing the average deviation of the twisted double plaquette from integer 
multiples of a, and then assign to Q the integer closest to Qt/ot. The factor a turns 
out to depend on (3, and approaches one with increasing /3; for example a 0.93 for 
(3 = 6 and a ~ 0.95 for (3 = 6.2, at the 20 th cooling step. This method eliminates 
the need for expensive, protracted cooling; usually fewer than 20 steps suffice (when 
using N(N — l)/2 subgroups) in order to observe a substantial convergence of the 
results, which appears to improve with increasing /3. The data that we will report 
are taken after 20 cooling steps. This cooling method to estimate Q is significantly 
less expensive than the overlap method, although less rigorous Q. 



To begin with, we report results from high-statistics MC simulations at 9 = 0. 
They required long MC runs, essentially because of the severe (apparently exponen- 



tial) critical slowing down of the topological modes | 15| , whose autocorrelation time 
(in units of sweeps) ranges from tq w 77 at (3 = 6 to tq « 540 at (3 = 6.2. In 
particular, over 40 million sweeps per value of /3 were produced on average for runs 
at 9 = 0. The results for the topological susceptibility Xi the coefficients 62 and 64, 
and the renormalization constant Z q , are reported in Table |l| The estimates of \ 
and 62 are consistent with earlier MC simulations, see, e.g., [0. The estimate of 
high-order coefficients, such as 64, turns out to be very hard in 9 = MC simulations, 
requiring huge statistics, due to the large cancellations in their expressions in terms 
of correlators at 9 = 0, see Eq (|9|), giving rise to relatively large statistical errors. 
The results for 6 4 are consistent with zero, suggesting the bound |6 4 | < 0.005, which 
will be improved below. The estimates of Z q may be compared with the one-loop 
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9 10 11 



Figure 2: Distribution of the ratio Qt/ot, for /3 = 6.2 configurations at 0l = and 9l 
(corresponding to 0j ps 1.5), after 20 cooling steps. 



8.7 



calculation |22| Z g = 1 — 5.448/3 + 0(/3 ), which underestimates it, and the es- 
timates by the heating method at p = 6: Z q = 0.18(1) and Z q = 0.16(1) 0, 
which are slightly larger than the values that we obtain in our MC simulations using 
Eq. (|1]). 

The MC simulations at 9l ^ are quite slower, due to the more complex struc- 
ture of the action terms involving a single link. The "force" corresponding to a given 
link now consists of additional contributions, beyond the standard sum of staples, 
as shown in Fig. |T[ Thus, the link updating procedure is slower by approximately a 
factor of three with respect to the 9l = case. In runs with 9l ^ 0, an average of 3 
million sweeps were produced for each value of (3 and 9l- 

The cooling algorithm is implemented in the same way as for the 6l = simu- 
lations. The presence of a 6l term with 6l > (9l < 0) leads to a preponderance 
of configurations with Q > (Q < 0). Fig. ^] shows the distributions of the ratio 
Qt/ot- of (3 = 6.2 cooled configurations at 6 L = and 0l = 8.7 (corresponding to 
6i = Z q 6 L « 1.5), after 20 cooling steps; we recall that the topological estimator Q 
is obtained by taking the integer closest to Qt/ot, a(f3 = 6.2) = 0.95. We note that 
these distributions cluster around integer values, also for rather large values of Qt/ot, 
both for 6 L = and 9l = 8.7, although the optimal value a = 0.95 is kept fixed. 
The autocorrelation time of the topological charge in the MC simulations at 9l ^ 
does not significantly change in units of sweeps with respect to that at 6 L = and 
the same (3. 

Figs. and | show MC results for the ratios (Q)e L / (Q 2 ) e=0 and (Q)e L / (Q 2 )e L , 
respectively, versus 0j = Z„9l- The MC data at different (3 values appear to follow the 
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Figure 3: The ratio {Q)e L / (Q 2 )e=o y s #i = Z q 9i,. Here scaling corrections are only visible 
for 0i > 2. The dashed line shows the curve Q% — 2b 2 8f with b 2 = -0.026, while the dotted 
line shows the simple linear behavior. 



same curve, providing evidence of scaling. Scaling corrections, which are expected 
to be 0(a 2 ), appear quite small, and tend to increase with increasing 0£. This good 
scaling behavior thus corroborates the existence of a nontrivial continuum limit for 
any value of 

In Table || we present the results of fits of the data for the ratio (Q)e L / (Q 2 )e=o, 
to the polynomial Ansatz 



{Q)e L 
(Q 2 ) 



0=0 



Z q 9 L - 2b 2 {Z q 9 L f + 3b 4 (Z q e L ) 5 + ... (24) 



The results turn out to be quite stable with respect to the maximum value of \6i\ 
allowed in the fits. They are in good agreement with the 6 = results presented in 
Table improving them significantly. In particular, they provide a much smaller 
bound on the value of 6 4 . The expected 0(a 2 ) scaling corrections in the estimates of 
b 2 i are taken into account by extrapolating the results at finite f3 fitting them to 

b 2i + ca((3) (25) 

where a is the string tension in unit of the lattice spacing (we use the data of 



Ref. P^[ ). This provides an accurate estimate of b 2 , that is 

b 2 = -0.026(3), (26) 

whose error includes statistical and systematic errors related to the extrapolation 
to the continuum limit and to the small differences of the fit results reported in 
Table [| This estimate of b 2 is clearly more precise than the one we may obtain 



-7- 



A 
V 



A 

v 



2.0 



1.5 



1.0 



0.5 



0.0% 



o p=5.9 
x p=6.0 
□ (3=6-2 



® x 



Figure 4: The ratio {Q)o L /{Q% vs { = Z q 8 L 



by a continuum-limit extrapolation of the data at 9 = reported in Table [l] using 
Eq. (p5|), which gives 6 2 = —0.029(7). Moreover, the results of Table || lead to a very- 
small bound for 64: 



|6 4 | < 0.001, 



(27) 



which is much smaller than the one obtained from the 9 = simulations. These 
results confirm that the coefficients b 2i of the expansion (|11]) are rapidly decreasing 
with increasing i, as expected by the large- N prediction b 2 i = 0(N~ 2% ) [[|. 

The results (p6|) and ( [27|) significantly improve earlier results for the expansion 
of the ground-state energy around 9 = 0, as obtained by 6* = MC simulations using 
different methods to estimate Q: b 2 = -0.023(7) from Ref. |]|, b 2 = -0.024(6) 
from Ref. ||17|| , and b 2 = —0.025(9) from Ref. |l8| , while no estimates for 64 had been 
reported. 

In conclusion, we have investigated the scaling behavior of the 4D SU(3) gauge 
theory in the presence of an imaginary 9 term, by MC simulations computing the 
topological properties at imaginary 9, i.e. 9 = —i9,i with real The numerical 
results for the topological charge provide a good evidence of scaling in the region 
\9i\ < 7r which we consider. The imaginary 9 dependence of the ground-state energy 
turns out to be well described by the first few nontrivial terms of the expansion 
around 9 = 0. Fits to polynomials provide a quite accurate estimate of b 2 and a 
very small bound on 6 4 , see Eqs. (p6|) and (|27D, which support the expected large- N 
scenario predicting b 2 i = 0(N~ 21 ). This study may be straightforwardly extended to 
other observables, to determine their 9 dependence, see Ref. ||32|| . Finally, we mention 
that, besides allowing more precise determinations of the 9 = expansion coefficients 
of the ground-state energy and other observables, the use of imaginary 9 values might 
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2.0 


0.1717(13) 


-0.027(2) 








2.7 


0.1717(12) 


-0.0265(13) 






0(9 5 L ) 


1.0 


0.1724(15) 


-0.027(4) 


-0.002(6) 






1.5 


0.1722(13) 


-0.027(4) 


-0.0003(14) 






2.0 


0.1722(13) 


-0.027(3) 


-0.0002(6) 






2.7 


0.1721(12) 


-0.026(3) 


-0.0000(2) 



Table 2: Results of some fits of the data of (Q)e L / (Q 2 )e=o to O(0\) odd polynomials. 
The fits use also the 9 = results of Table |l] and consider only data up to a maximum 
value |6>i|max of \6i\, leading to an acceptable x 2 i X 2 /d°f < 2. 



turn out useful in the effort to overcome the problem of the dramatic critical slowing 



down of the topological modes |4j, 0, 0, [34], f35|, by performing parallel tempering 
simulations j36j with a set of imaginary 9 values including 8 = 0, which provides an 
exact MC algorithm for the model at 9 = 0. This approach is largely used in the 
MC simulations of spin-glass models which are affected by analogous critical slowing 
down problems. 
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